3D plotting



In [1]:

    
from mayavi import mlab









    



WARNING:traits.has_traits:DEPRECATED: traits.has_traits.wrapped_class, 'the 'implements' class advisor has been deprecated. Use the 'provides' class decorator.



In [2]:

    
dir(mlab)









    Out[2]:





['ETSConfig',
 '__builtins__',
 '__doc__',
 '__file__',
 '__name__',
 '__package__',
 'animate',
 'axes',
 'barchart',
 'clf',
 'close',
 'colorbar',
 'contour3d',
 'contour_surf',
 'draw',
 'figure',
 'flow',
 'gcf',
 'get_engine',
 'imshow',
 'mesh',
 'move',
 'options',
 'orientation_axes',
 'outline',
 'pipeline',
 'pitch',
 'plot3d',
 'points3d',
 'quiver3d',
 'roll',
 'savefig',
 'scalarbar',
 'screenshot',
 'set_engine',
 'show',
 'show_engine',
 'show_pipeline',
 'start_recording',
 'stop_recording',
 'surf',
 'sync_camera',
 'sys',
 'test_barchart',
 'test_contour3d',
 'test_contour3d_anim',
 'test_contour_surf',
 'test_fancy_mesh',
 'test_flow',
 'test_flow_anim',
 'test_imshow',
 'test_mesh',
 'test_mesh_sphere',
 'test_mesh_sphere_anim',
 'test_molecule',
 'test_plot3d',
 'test_plot3d_anim',
 'test_points3d',
 'test_points3d_anim',
 'test_quiver3d',
 'test_quiver3d_2d_data',
 'test_simple_surf',
 'test_simple_surf_anim',
 'test_surf',
 'test_triangular_mesh',
 'text',
 'text3d',
 'title',
 'triangular_mesh',
 'vectorbar',
 'view',
 'wxversion',
 'xlabel',
 'yaw',
 'ylabel',
 'zlabel']



In [3]:

    
mlab.test_mesh()
mlab.show()



In [4]:

    
mlab.test_triangular_mesh()
mlab.show()



In [5]:

    
mlab.test_quiver3d()
mlab.show()



In [6]:

    
import numpy as np

x = np.linspace(0, 2*np.pi, 128)

f = (np.sin(  x[:, np.newaxis, np.newaxis]+3*x[np.newaxis, :, np.newaxis]) +
     np.cos(2*x[np.newaxis, np.newaxis, :]+2*x[np.newaxis, :, np.newaxis]))

mlab.contour3d(f.T)
mlab.show()



In [7]:

    
f = (np.sin(  x[:, np.newaxis, np.newaxis]+3*x[np.newaxis, :, np.newaxis]) +
     np.cos(2*x[np.newaxis, np.newaxis, :]+2*x[np.newaxis, :, np.newaxis]))

g = (np.sin(2*x[:, np.newaxis, np.newaxis]+3*x[np.newaxis, :, np.newaxis]) +
     np.cos(  x[np.newaxis, np.newaxis, :]+2*x[np.newaxis, :, np.newaxis]))

h = (np.sin(3*x[:, np.newaxis, np.newaxis]+  x[np.newaxis, :, np.newaxis]) +
     np.cos(  x[np.newaxis, np.newaxis, :]+2*x[np.newaxis, :, np.newaxis]))

#help(mlab.flow)
mlab.flow(f, g, h)
mlab.show()



In [10]:

    
def generate_data_3D(
        n,
        p = 1.5,
        seed = None):
    assert(n % 2 == 0)
    np.random.seed(seed)
    a = np.zeros((n, n, n//2+1), dtype = np.complex)
    a[:] = np.random.randn(*a.shape) + 1j*np.random.randn(*a.shape)
    k, j, i = np.mgrid[-n//2+1:n//2+1, -n//2+1:n//2+1, 0:n//2+1]
    k = np.roll(k, n//2+1, axis = 0)
    j = np.roll(j, n//2+1, axis = 1)
    kk = (k**2 + j**2 + i**2)**.5
    a /= kk**p
    a[0, :] = 0
    a[:, 0] = 0
    ii = np.where(kk == 0)
    a[ii] = 0
    ii = np.where(kk > n/3.)
    a[ii] = 0
    return np.fft.irfftn(a)

f = generate_data_3D(32, p = 2)

mlab.contour3d(f.T)
mlab.show()



In [14]:

    
import itertools
points = itertools.product([0, 1], repeat = 3)
for p in points:
    print(p)
    mlab.plot3d(np.array([0, p[0]]),
                np.array([0, p[1]]),
                np.array([0, p[2]]),
                color = (1, 0, 0))
mlab.show()









    



(0, 0, 0)
(0, 0, 1)
(0, 1, 0)
(0, 1, 1)
(1, 0, 0)
(1, 0, 1)
(1, 1, 0)
(1, 1, 1)

Lorenz equations

\begin{align} \frac{d}{dt} x &= \sigma (y - x) \\ \frac{d}{dt} y &= x (\rho - z) - y \\ \frac{d}{dt} z &= x y - \beta z \end{align}

chaotic for \begin{equation} \sigma = 10, \beta = 8/3, \rho = 28 \end{equation}



In [33]:

    
sigma = 10.
beta  = 8./3
rho   = 28.

def lorenz_rhs(x):
    return np.array([sigma * (x[1] - x[0]),
                     x[0]*(rho - x[2]) - x[1],
                     x[0]*x[1] - beta*x[2]])

def Euler(rhs, dt = 0.5, N = 8, x0 = None):
    x = np.zeros((N+1,) + x0.shape,
                 dtype = x0.dtype)
    x[0] = x0
    for t in range(N):
        x[t+1] = x[t] + dt*rhs(x[t])
    return x

def cRK(rhs, dt = 0.5, N = 8, x0 = None):
    x = np.zeros((N+1,) + x0.shape,
                 dtype = x0.dtype)
    x[0] = x0
    for t in range(N):
        k1 = rhs(x[t])
        k2 = rhs(x[t] + 0.5*dt*k1)
        k3 = rhs(x[t] + 0.5*dt*k2)
        k4 = rhs(x[t] + dt*k3)
        x[t+1] = x[t] + dt*(k1 + 2*(k2 + k3) + k4)/6
    return x

lsol = cRK(lorenz_rhs,
             dt = 2.**(-8),
             N = 1024,
             x0 = 10+np.random.random((3, 32)))

for traj in range(32):
    mlab.plot3d(lsol[:, 0, traj],
                lsol[:, 1, traj],
                lsol[:, 2, traj],
                tube_radius = None,
                color = (traj*1./32, 0, 1 - traj*1./32))
mlab.show()



In [ ]: